/***********************************************************************
 * This program implements a simple molecular dynamics simulation,
 *   using the velocity Verlet time integration scheme. The particles
 *   interact with a central pair potential.
 *
 * Author:   Bill Magro, Kuck and Associates, Inc. (KAI), 1998
 *
 * Parallelism is implemented via OpenMP directives.
 ***********************************************************************/

#define SEQ_TIME 0

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <sys/types.h>
#include <sys/time.h>
#include <unistd.h>

#ifndef RAND_MAX
#define RAND_MAX 0x7fff
#endif

#define ndim 3
#define nparts 8192
#define nsteps 10

typedef double real8;

typedef double vnd_t[ndim] ;

double t() {
    struct timeval tv;
    gettimeofday(&tv, ((void *)0));
    return (double)tv.tv_sec + (double)tv.tv_usec/1000000.0;
}

/* statement function for the pair potential and its derivative
   This potential is a harmonic well which smoothly saturates to a
   maximum value at PI/2.  */

real8 v(real8 x) {
  if (x < M_PI_2) 
    return pow(sin(x), 2.0);
  else
    return 1.0;
}

real8 dv(real8 x) {
  if (x < M_PI_2) 
    return 2.0 * sin(x) * cos(x);
  else
    return 0.0;
}


/***********************************************************************
 * Compute the forces and energies, given positions, masses,
 * and velocities
 ***********************************************************************/
void compute(vnd_t *f)
{
  int j, k;
  vnd_t rij;
  
  /* The computation of forces and energies is fully parallel. */
#pragma omp target device(mpi)
#pragma omp parallel for default(shared) private(rij)
    /* compute potential energy and forces */
    for (j = 0; j < 9; j++) {
	    for (k = 0; k < 9; k++) {
		    f[j][k] = rij[k];
		}
      }
}
       
/******************
 * main program 
 ******************/
int main (int argc, char **argv) {
  
  /* simulation parameters */
  vnd_t force[nparts];
  int i;
  
  for (i = 0; i < ndim; i++)
    box[i] = 10.0;
    
  /* compute the forces and energies */
  compute(force);

  return (0);
}

